Objeto de patrón de puntos

Ahora podemos usar la función ppp, en spatstat, para crear el objeto de patrón de puntos:

acc.ppp <- ppp(x=st_coordinates(act.sp.401)[,1],
               y=st_coordinates(act.sp.401)[,2],
               window=window)

Intensidad y densidad

Una información crucial que necesitamos cuando tratamos con patrones de puntos es una definición cuantitativa de la distribución espacial, es decir, cuántos eventos tenemos en una ventana predefinida. El índice para definir esto es la Intensidad, que es el número promedio de eventos por unidad de área.

Podemos calcular la intensidad de la siguiente manera:

acc.ppp$n/sum(st_area(heredia))

0.00005973704 [1/m^2]

Dato de intensidad

El numerador es el número de punto en el objeto ppp; mientras que el denominador es la suma de las áreas de todos los polígonos. En Heredia, la intensidad promedio de accidentes es de 59.737 por kilómetro cuadrado.

Observación

La intensidad puede ser constante a lo largo de la ventana de estudio, en ese caso, en cada metro cuadrado encontraríamos el mismo número de puntos, y el proceso sería uniforme o homogéneo. La mayoría de las veces, la intensidad no es constante y varía espacialmente a lo largo de la ventana del estudio, en ese caso el proceso no es homogéneo. Para procesos no homogéneos necesitamos una forma de determinar la cantidad de variación espacial de la intensidad.

Conteo de cuadrantes

Hay varias formas de tratar este problema, un ejemplo es el conteo de cuadrantes, donde el área se divide en rectángulos y se cuenta el número de eventos en cada uno de ellos:

plot(acc.ppp,pch="+",cex=0.5,
     main="Heredia: Accidentes de tránsito")
plot(quadratcount(acc.ppp, nx = 4, ny = 4), add=T,col="blue")

Accidentes por distrito

Esta función es buena para ciertos conjuntos de datos, pero en este caso realmente no tiene sentido usar el recuento de cuadrantes, ya que las áreas que crea no tienen ningún significado en la realidad. Por ejemplo, sería mucho más valioso extraer el número de accidentes por distrito. Para hacer esto necesitamos usar un bucle e iterar a través de los polígonos:

Local.Intensidad <- data.frame(heredia=factor(),Number=numeric())
for(i in unique(heredia$distrito)){
  sub.pol <- heredia[heredia$distrito==i,]
  z=data.frame(cbind(x=acc.ppp$x, y= acc.ppp$y))
  coordinates(z)<-c("x","y")
  sub.ppp <- ppp(x=z$x, y=z$y,window=as.owin(sub.pol))
  Local.Intensidad <- rbind(Local.Intensidad,
                            data.frame(heredia=factor(i,
                                                      levels=heredia$distrito),
                                       Number=sub.ppp$n))
}

Frecuencia

Se establecen las frecuencias de los accidentes por distrito.

summary(Local.Intensidad)

          heredia      Number     
 Ulloa        :1   Min.   :152.0  
 San Francisco:1   1st Qu.:275.8  
 Heredia      :1   Median :368.5  
 Mercedes     :1   Mean   :372.5  
                   3rd Qu.:465.2  
                   Max.   :601.0  

sum(Local.Intensidad$Number)

[1] 1490

Código de Gráfico de frecuencia

library(plotrix)
colorScale <- color.scale(Local.Intensidad[order(Local.Intensidad[,2]),2],
                          color.spec="rgb",
                          extremes=c("green","red"),alpha=0.8)

df_sorted <- Local.Intensidad %>% arrange(desc(Number))

p<-ggplot(data=df_sorted, aes(x=heredia, y=Number)) +
  geom_col(aes(color = heredia), fill = "white")+
  scale_fill_manual(colorScale) +
  geom_text(aes(label = Number), vjust = -0.3) +
  theme_minimal()

Gráfico de frecuencia

p

Figura 2: Heredia: Frecuencia de accidentes de tránsito por distrito.

Densidad de Kernel

Otra forma en la que podemos determinar la distribución espacial de la intensidad es mediante el uso del suavizado del kernel (Diggle 1990; Berman 1986; Berman, Scientific, y Diggle 1989; Bivand, Pebesma, y Gómez-Rubio 2013). Dicho método calcula la intensidad de forma continua en toda el área de estudio. Para realizar este análisis en R necesitamos definir el ancho de banda de la estimación de densidad, que básicamente determina el área de influencia de la estimación. No hay una regla general para determinar el ancho de banda correcto.

Densidad de Kernel

En general, si h es demasiado pequeño, la estimación es demasiado ruidosa, mientras que si h es demasiado alta, la estimación puede pasar por alto elementos cruciales del patrón de puntos debido a exceso de suavisado (Scott 2009). En spatstat, las funciones bw.diggle, bw.ppl y bw.scott se pueden usar para estimar el ancho de banda según los métodos de diferencia. Podemos probar cómo funcionan con nuestro conjunto de datos utilizando el siguiente código:

Densidad de Kernel 01

plot(density.ppp(acc.ppp, 
                 sigma = bw.diggle(acc.ppp),edge=T),
     main=paste("h =", round(bw.diggle(acc.ppp),2)))

Figura 3: Heredia: Densidad de accidentes de tránsito por distrito.

Densidad de Kernel 02

plot(density.ppp(acc.ppp, 
                 sigma = bw.ppl(acc.ppp), edge=T),
     main=paste("h =", round(bw.ppl(acc.ppp),2)))

Figura 4: Heredia: Densidad de accidentes de tránsito por distrito.

Densidad de Kernel 03

plot(density.ppp(acc.ppp, 
                 sigma = bw.scott(acc.ppp)[2], edge=T),
     main=paste("h =", round(bw.scott(acc.ppp)[2],2)))

Figura 5: Heredia: Densidad de accidentes de tránsito por distrito.

Densidad de Kernel 04

plot(density.ppp(acc.ppp, 
                 sigma = bw.scott(acc.ppp)[1],edge=T),
     main=paste("h =", round(bw.scott(acc.ppp)[1],2)))

Figura 6: Heredia: Densidad de accidentes de tránsito por distrito.

Muchas gracias